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SUMMARY 


Further computational experiments have been conducted to study the char- 
acteristics of flow reversal and separation in unsteady boundary layers. 
One set of calculations was performed using the first-order, time- 
dependent turbulent boundary-layer equations, and extended earlier work 
by Nash and Patel to a wider range of flows. Another set of calcula- 
tions was performed for laminar flow using the time-dependent Navier- 
Stokes equations. 

The results of the calculations confirm previous conclusions concerning 
the existence of a regime of unseparated flow, containing an embedded 
region of reversal, which is accessible to first-order boundary-layer 
theory. However, certain doubts are cast on the precise nature of* the 
events which accompany the eventual breakdown of the theory due to 
singularity onset. The earlier view that the singularity appears as the 
final event in a sequence involving rapid thickening of the boundary 
layer and the formation of a localized region of steep gradients, is 
called into question by the present results. It appears, first, that 
singularity onset is not necessarily preceded by rapid boundary-layer 
thickening, or even necessarily produces immediate thickening. Further- 
more, the formation of a region of steep gradients could not be repro- 
duced in the solutions of the Nav i er-Stokes equations, and may, itself, 
prove to be a feature of first-order bounda ry- 1 ayer theory and not part 
of a more complete description of the flow. 
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INTRODUCTION 


Extensive calculations have been performed, over the last several years, 
to study the characteristics of time-dependent turbulent boundary layers. 
These calculations have taken the form of numerical experiments, yielding 
data which -- although inferior in obvious respects to reliable wind-tunnel 
results — have nevertheless served to fill the void which still exists 
because of the scarcity of comprehensive measurements. Indeed, much of 
what is known at this time about unsteady boundary layers has been derived 
from theory rather than from measurement, and, while the need for good 
wind-tunnel experiments to confirm the calculated results (or at least to 
validate the theoretical models) remains as acute as ever it was, the 
information content of calculations seems likely to "stay ahead" of that 
of measurements for some time to come. 

The most recent work on unsteady turbulent boundary layers, done by the 
present author and his co-workers, has been concerned with the effects 
of time-dependence on reversal and separation onset. It was shown in 
References [1,2] that time-dependence results in a delay of reversal onset, 
and that this delay could not be related in any simple manner to an alle- 
viation of the pressure gradients; indeed, some measure of delay was 
observed even when the pressure gradients were augmented by the effects 
of the unsteadiness. It was shown further, that the reversal point in an 
unsteady turbulent boundary layer is not a singular point -- an observa- 
tion which confirmed earlier statements which had been made concerning 
laminar flow [3, ^,5,6] -- and that the boundary layer remains thin even 
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though reversal had taken place [1], Subsequent work was directed to 
studying a class of unsteady turbulent boundary layers In which a region 
of embedded reversal was allowed to develop but in which, under certain 
conditions, a singularity occurred sometime later [7,8]. This singu- 
larity is related to the final separation of the boundary layer, in the 
sense of detachment of the outer flow from the body surface, and repre- 
sents the limit of validity of first-order boundary- 1 ayer theory. However, 
the regime between reversal and separation appears to be accessible to 
first-order theory [7], and the calculations which were performed in this 
regime are believed to be meaningful. 

The present work, which forms a sequel to that of Reference [7] , has two 
main objectives: 

1. to extend the calculations of unsteady turbulent boundary layers 
to a wider range of flow conf igurat ions , and 

2. to try to elucidate the mechanism of separation onset by perform- 
ing higher-order calculations which would not break down as the 
result of the development of a singu lari ty. 

In pursuit of the first objective, calculations have been performed for 
two additional types of flow: flows in which an unsteady retardation is 

first imposed, but where the external velocity distribution is subse- 
quently frozen allowing a relaxation towards steady-state conditions, and 
oscillatory flows in which reversal occurs during part of the cycle. The 
results of these calculations are presented in Part I of this report. 

The second objective was pursued by programming a simple method for solv- 
ing the unsteady Nav i er-S tokes equations, and performing some calculations 
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for laminar flow. It was recognized that the characteristics of turbu- 
lent flows would not be represented appropriately by this means, but it 
was hoped that the mechanism of separation would have sufficient generality 
to make these laminar calculations useful. The results of this study are 
presented in Part II of the report. 
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PART I. TURBULENT FLOW CALCULATIONS 


1.1 NATURE OF THE FLOWS CONSIDERED 

The turbulent boundary-layer studies reported here were carried out for 
incompressible, time-dependent flow over a two-dimensional surface of 
large or infinite radius. As in the earlier work [1,7], orthogonal 
coordinates are erected on the surface, with y measured normal to it and 
x measured along the surface from some origin where the boundary layer is 
already turbulent and of known properties: correspond i ng to steady, con- 

stant-pressure conditions. The main features of the calculation method 
are reviewed in Appendix A. The calculations relate to the boundary layer 
developing over a “plate" of chord c which extends downstream from the 
origin. 

The external velocity: U^, over the plate, is assumed to vary in a pre- 

scribed manner with x and time, t. Specifically, it is assumed that 

U = U , for t < 0 and all x (I) 

e o — 

U 

= 1 - {1 - f (t) }, for t >_ 0 and 0 < x < x q (2) 

o o 

U x - x 

- 1 - — {l - f(t)}, for t > 0 and < x <_ c (3) 

o o 

where U is some reference velocity, and x and x are prescribed values 
o o l 

of x: x = 5c/7 ( = 0 . 71 ^ 0 ), X = 2x . In Equations (2,3), f(t), which 

o I o 

is the value of U /U at x = x , is chosen to be a function of time, with 




the restriction that f = 1 at t = 0 in order to satisfy Equation (1). 

The flows discussed in Reference [7] corresponded to 

f = 1 - cox t/c, (k) 

o 

where to = constant, which imposed a distortion of the external flow which 
was monotonic in time. 

Here, the emphasis is on patterns of distortion of the external velocity 
field which are of two other types: 

(1) patterns in which an initial distortion, extending over a finite 

' interval (0 t _< t^) , is followed by freezing of the external 

velocity field, 

(2) oscillatory external velocity fields. 

The "frozen 11 flows were generated by taking 

f = 1 - (1 - f f )t/t f , for 0 < t < t f (5) 

f = f f , for t t f (6) 

Values of f^ = 0 and 0.5 were used, together with values of t ^ corre- 
sponding to U^t^/c = 0.2, 1.0, and 2.0. 

The oscillatory flows were of a triangular waveform, and were generated 
by taking 

t - t 

f = 0.5 - A 4{— 2. - 1} (7) 

P 

over the first half of the cycle: t < t < t + t /2, and 

9 A — — n n 
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f = 0.5 + A{4 


( 8 ) 



over the second half of the cycle: t +t/2<t<t. Here, t is the 

o p — — o p 

period of the oscillation, and t Q refers to the beginning of the cycle. 

It will be noted that the mean value of f is 0.5, except* for the first 
half-cycle of the motion (0 <_ t <_ t p /2) where an initial transient sat- 
isfies the requirement: f = 1 for t = 0; this transient is defined by 

f = 1 - (1 + 2A) . (9) 

P 

The period of the motion, which in dimensionless terms is U t /c, was 
r ’ op 

varied between 2.0 and 8.0 in the calculations. The constant A, in Equa- 
tions (7,8,9), controls the amplitude of the oscillatory external velocity 
distribution; values between 0 and 0.5 were chosen for the calculations. 
When A = 0, the flow degenerates from an oscillatory to a "frozen" flow 
(as defined above) with f^ = 0.5 and t^ = t^/2. 

When A = 0.5, in the oscillatory calculations, the external flow stagnates 
momentarily (i.e. = 0), at the point on the surface: x = x q , once 

during every cycle: namely, at t = t Q + t p /2. = in the 

"frozen flows," stagnation of the outer flow occurs at x = x q for all 
times subsequent to t = tp Stagnation does not occur at al 1 if A < 0.5, 
in the oscillatory flows, and i f f ^ > 0 in the frozen flows. 

Figure 1 illustrates the features of the two types of flow, and provides 
a graphical definition of the parameters fp tp t^ and A. 

A 

If A = 0.5 the mean value of f is 0.5 even during the first half-cycle. 
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1.2 COMPUTATIONAL EXPERIMENTS 


Frozen Flows 

A number of calculations were done for flows in which the external 
velocity distribution was frozen, at t = t^, following an initial dis- 
tortion starting at t = 0. 

Figures 2 through 7 show the spacial distributions of wall shear stress 
and displacement thickness, respectively, at various time levels, for 
the following cases: 

U t./c 

o f f 

0.2 0 

2.0 0 

2.0 0 . 5 . 

In the first two cases, a region of reversed flow, of substantial extent, 
has already formed on the plate by the time (t = t^_) the external flow 
is frozen. In the third, reversal has not yet taken place at the time 
of freezing but occurs later. With f^. = 0, stagnation of the outer flow 
occurs at x = x , and t > t r , and the displacement thickness takes on a 
locally infinite value. However, this infinity reflects the division by 
zero in the definition of displacement thickness and does not necessarily 
indicate the presence of a singularity in the solution of the boundary- 
layer equations. 

Various criteria were described in Reference [7] for identifying the onset 
of a singularity: steepening of the gradients of wall shear stress and 

displacement thickness, versus x, the development of locally high dis- 
placement-thickness maxima, and breakdown of the boundary- 1 ayer approxima- 

JL 

tioned. A surrogate criterion was also suggested: 6“ = 0.1c, for flows 
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of the present type and at these Reynolds numbers, which appeared to 
correlate with the other criteria, and which helped to pinpoint the loca- 
tion of the locally high gradients between node points in the calcula- 
tion. This surrogate criterion will be used here also, for the moment, 
although its usefulness is called into question later in the discussion. 

In any case, it is necessary to repeat the cautionary statements made in 
Reference [7] to the effect that the real measure of singularity onset is 
steepening of the velocity gradients rather than the attainment of some 
arbitrary value of 6*, and that the surrogate only has validity insofar 
as it helps to locate the formation of the steep gradients. 

Figures 2 through 5, which relate to the two cases where f ^ = 0, indicate 
that incipient singular conditions form at a time corresponding to U Q t/c=2.0 
for L^t^/c = 0.2, and at a time corresponding to U Q t/c = 3-0 for l^t^/c = 2.0. 
Thus the lower initial rate of distortion of the outer flow, in the latter 
case, results in a delay in the final onset of the singularity. Figures 8, 

9, in which the wall shear stress and the displacement thickness are plotted 
versus time, emphasize this delay. The small dip in one of the curves of 
6* versus time (Figure 9) results from a forward movement of the 6* maximum 
on the plate. It is interesting to note that, for these two cases with f^=0, 
there is a measure of similarity in the development of each solution towards 
singular conditions: the two sets of results (for or 6 ) would roughly 

coincide if they were plotted against t-t s> where t^ refers to the onset of 

the singularity. This similarity does not extend to the case where f ^ = 0.5 
(Figures 6,7), where the approach to the singularity is somewhat slower 
(Figures 8,9). The slower approach to singular conditions presumably has to 
do with the smaller gradient of with respect to x, when f^ 55 0.5, but direct 

comparison of the results for this case with those for f ^ = 0 is difficult 
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because the singularity forms further aft when f^ = 0.5: at around 

x = 0.6c, rather than x = 0.25c with f^. = 0. 

In Reference [7] a semi -empi r i ca 1 model was proposed for predicting the 
onset of the singularity. The procedure was to compare the velocity of 
penetration of the reversal point, upstream into the oncoming boundary 
layer, with an average "convection" velocity in the reversed-f low region. 
For want of a better definition, this convection velocity was taken to 
be one-half of the maximum negative velocity in the region at the particu- 
lar time level. The prediction criterion proved useful in correlating 
the results of Reference [7], and it was therefore of interest to see 
whether it would be valid here. To this end, estimates were made of the 
penetration velocity, for the three frozen flows considered here, by plot- 
ting the movement of the point of flow reversal (Figure 10), and differ- 
entiating the curves graphically. The values obtained in this manner 
are shown in Figure 11, and compared with data for the convection velocity, 
as defined above. 


The points of intersection of the curves corresponding to the pairs of 
velocities yields the predicted times at which the singularity is sup- 
posed to form. These predictions are compared with the approximate 
observed values (indicated by 6" = 0.1c) in the following table: 


Case: 

U t,/c f. 

of f 


Predicted onset Observed onset 


0.2 

2.0 

2.0 


0 

0 

0.5 


U t/c 

0 

1 .2 
2.2 
3.9 


U t/c 
o 

2.0 

3.0 

4.8. 


9 



It will be seen that there is a constant discrepancy in U Q t/c, of about 
0.8, between the predictions and the observations, and the criterion has 
no more than qualitative validity for these results. The reason why it 
does not apply here is not clear, although some suggestions can be made. 

It might be suggested that the assumed ratio of convection velocity to 
maximum negative velocity needs to be reduced. However the d i screpanc i es 
are too large to be resolved in such a manner: the ratio would have to be 
reduced virtually to zero to make the predictions and observations agree. 
Indeed, for the present frozen flows, the singularity does not seem to 
form until the forward movement of the reversal point has ceased: i .e. 

until the penetration velocity has fallen to zero. In the flows discussed 
in Reference [ 7 ], the singularity was observed to develop while the 
reversal point was still moving forward. Another suggestion, antici- 
pating some of the later discussion in this report, is that the real 
singularity, in the sense of steepening of the gradients, does form at 
the predicted times, and that the discrepancies simply measure the failure 
of the surrogate criterion to pinpoint its formation. Since the surro- 
gate criterion was used in the derivation of the heuristic model, in 
Reference [ 7 ], the implication must be that the rates of increase of 6“ 
with time were so much greater than they are here that any differences 
in timing between true singularity onset and the condition: <5'' = 0.1c, 
were masked. It is true that the rates of increase of 6* were generally 
higher, in Reference [ 7 ], than they are in the present flows, and so the 
suggestion has some merit. On the other hand the discrepancies, between 
the predicted and the observed singularity onset are rather larger than 
can be explained entirely in this way, and for the moment the conclusion 
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must remain that the heuristic model, proposed in Reference [7], is not 
sufficiently general to make quantitative predictions outside the class 
of flows for which it was originally developed. 
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Oscillatory Flows 


(a) Two sets of calculations were performed for oscillatory flows; the 

first set was concerned with the effects of variation of frequency, and 

included the following cases: 

U t /c A 

o P 

2.0 0.5 

1*.0 0.5 

8.0 0.5- 

The function f(t), corresponding to each of these cases, is plotted in 
Figure 12. The pressure gradient over the plate: 3p/3x, is given by 


I lE 

p 9x 



31 ) 


U 


e 3x 


( 10 ) 


and varies with both x and t. The pressure gradient has a discontinuity 
at x = x q ; however, the value at x = x q is representative of conditions 
over the forward portion of the plate: 0 <_ x _< x q , and can be expressed 

as a function of f and its derivative f 1 (= df/dt) by: 


Vzlt* f(l ' f) ' f ' ■ <"> 

pU o o 

o 

Now, since from Equations (7-8), f' = -2/t p , over the first half of the 
cycle, and - 2/t^, over the second half, 8p/8x at x = x q can be determined 
as a function of time, and is shown in Figure 13 for the three cases 
defined above. 


12 



The distributions of wall shear stress and displacement thickness, for 
the three cases, are presented in Figures 14 through 19- The emphasis, 
here, is on demonstrating the failure of the boundary layer to reach a 
stable oscillatory condition. By contrast, when reversal does not occur 
during part of the cycle, stable conditions are reached after a few cycles 
(see Reference [9]). For U^t^/c « 2.0 (Figures 14,15), the distributions 
at stages through the first cycle are compared with those at correspond- 
ing stages during the second cycle. It will be seen that the wall shear 
stress is typically lower, and the displacement thickness typically 
higher, during the second cycle than the first. The differences are most 
conspicuous over the range of x corresponding to the reversed-f 1 ow region. 
This trend continues through the third cycle, and is evidently associa- 
ted with the gradual approach towards eventual separation and the forma- 
tion of a singularity. Nevertheless, during the latter part of the first 
cycle and the early part of the second, the reversed-f low region shrinks 
and vanishes, and forward flow is temporarily reestablished over the 
whole of the plate. This intermittent behavior which will be discussed 
more fully later, is repeated for at least two more cycles. 

For the intermediate frequency: U^t^/c = 4.0 (Figures 16,17), reversed 

flow persists throughout the first cycle and a singularity forms towards 
the end of the cycle. It will be recalled from Figure 13 that 3p/3x < 0 
during the second half of the cycle, and so we have the interesting 
situation here in which reversed flow exists on the plate in spite of a 
favorable pressure gradient. With l^t^/c = 8.0, a singularity forms 
during the first half-cycle and forward flow is never recovered once 
reversal is established. 
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The movement of the reversal point, with time, is plotted in Figure 20, 

for the three cases. In all three cases, reversal first occurs near 

x = x q and then moves rapidly upstream with increasing time. At the 

highest frequency considered (U^t^/c = 2.0), the reversal point reaches 

a furthest upstream position at the mid point of the cycle and then moves 

aft again, with reversed-f low vanishing at a point on the plate somewhat 

upstream of x . At the two lower frequencies (U t /c = 4.0 and 8.0), the 
o o p 

reversal point never moves appreciably aft from its furthest-upstream 
position. The approximate times at which 6'' first reaches the value 0.1c 
(other than when 0) are shown on the curves in Figure 20. 


(n order to observe an intermediate behavior, between that for U t /c=2.0 

o p 

in which reversed flow temporarily vanishes during part of the cycle, and 

that for U t /c = 4.0 when it persists indefinitely, an additional calcu- 
o p r 7 

lation was made for U t /c = 2.8 and the results are shown by the dashed 

o p 

curve in Figure 20. Here, forward flow is reestablished for a short 
time, at the end of the first cycle, but reversed flow persists during 
the second cycle with 5* reaching the value 0.1c about half-way through 
the second cycle. 


Figures 21, 22 show the variation, with time, of the wall shear stress 
and displacement thickness at two neighboring positions on the plate, 
for this last case. The displacement thickness reaches 0.1c at about 
U^t/c = 4.0, although steep gradients of wall shear stress occur some- 
what earlier. It is of interest to note the pronounced phase shift 
between the displacement thickness and the external velocity gradients 
(which are symmetrical about the mid point of the cycle) even during 
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the first cycle when the flow is still far from developing a singularity. 
The magnitude of this phase lag is considerably greater than the lags 
found in Reference [9l for flows which did not suffer reversal during 
part of the cycle. 

The limited data presented here suggest that there is some connection 
between the persistence of reversal during the cycle, even when the 
pressure gradient becomes negative, and the appearance of the singular- 
ity. This may simply mean that a singularity can only develop when 
reversed flow is present; but, if so these results are at variance with 
suggestions which have been made that a singularity can occur upstream 
of reversal. Alternatively, the incipient formation of a singularity 
may somehow inhibit the shrinkage and disappearance of the reversed-f low 
region which would otherwise take place. 

(b) The second set of calculations for oscillatory flow was concerned 
primarily with the effects of variation of amplitude, and included the 
following cases: 


U t /c 
o p 

A 

2.0 

0 

2.0 

0.1 

2.0 

0.3 

2.0 

0.5 

4.0 

0 

4.0 

0.1 

4.0 

0.3 

4.0 

0.5. 


The two cases with A = 0.5 are the same as two of the oscillatory flows 

mentioned above, and the case: A = 0, U t /c - 4.0 is the same as the 

’ o p 

frozen flow, where U t r /c = 2.0 and f r = 0.5, which was discussed 

o t i 
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earlier. The function f(t) is plotted in Figure 23, for the cases where 

U t /c = 2.0, and in Figure 32 for the cases where U t /c = 4.0. The 
op op 

corresponding plots of the pressure gradient, at x = , are shown in 

Figures 24 and 33, respectively. 

The distributions of wall shear stress and displacement thickness, over 
the plate, for the cases where U^t^/c = 2.0, are shown in Figures 25, 

26 for a time corresponding to U Q t/c = 4.0: i.e., at the end of two com- 
plete cycles, and in Figures 27, 28 for U^t/c = 5.0, which corresponds 
to the mid point of the third cycle. In Figures 29, 30 the wall shear 
stress and displacement thickness at one point on the plate, are plotted 
versus time, and Figure 31 shows the movement of the reversal point with 
time. 

It will be noted from Figure 29 that, for A = 0 and 0.1, reversal does 
not occur until some time during the second cycle. However, once it 
does take place, it persists until a singularity develops in the third 
cycle (Figure 28). In contrast, for A = 0.3 and 0.5, reversal appears 
during the first cycle but subsequently vanishes again, later in the 
cycle, leaving forward flow reestablished over the whole plate. This 
intermittent formation and subsequent disappearance of reversed flow is 
repeated for at least two more cycles (Figure 29). The displacement 
thickness reaches a value of 0.1c at around the mid point of the third 
cycle, when A = 0.3, but there is no evidence of singular behavior any- 
where in the first three cycles when A = 0.5 (the calculation was not 
continued beyond U Q t/c = 6: the end of the third cycle); the large values 
of 6* at U Q t/c = 5 (Figure 28) result from the stagnation of the outer 
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flow and do not indicate the onset of a singularity in the solution. The 
difference in behavior between the two low-amplitude cases, in which 
reversal persists throughout the second and subsequent cycles, and the 
two high-amplitude cases, where reversal is intermittent, correlates 
with the sign of 3p/3x during the second half of the cycle (Figure 24): 
positive for A = 0, 0.1 and negative for A = 0.3, 0.5. On the other hand 
it has already been observed that reversal can persist even when the pres- 
sure gradient is favorable, and, clearly, the mechanism which controls the 
appearance and disappearance of the reversed-f low region is not simply a 
function of the pressure gradient. Indeed, as already mentioned, it was 
shown in Reference [1] that time-dependence caused a delay in the first 
appearance of reversal even when 3U/3t made a positive contribution to 
the pressure gradient. 

Figures 34 through 36 present data for the lower frequency calculations 

(U t /c = 4.0). Here, it will be seen that forward flow is not recovered 
o p 

for any of the amplitudes considered; once formed, reversal persists 
until 6* reaches 0.1c near the end of the first cycle. Specifically, 
singularity onset (according to this criterion) occurs earliest for the 
high-amplitude case: at U Q t/c = 3-7, and latest for A = 0: at about 

U Q t/c = 4.7. This order of onset, in terms of the values of A, is oppo- 
site to the order observed at the higher frequency. The difference is 
probably associated with the failure of the reversal point to move aft, 
during the second half of the cycle, when U Q t^/c = 4.0 (Figure 36). At 
the lower frequency, the reversal point moves forward, during the first 
half-cycle, and subsequently remains near its furthest- forward position, 
leaving an extensive region of reversed flow, for the duration of the 
motion. It is interesting to note the wide variation in the terminal 
positions of the reversal point, over the range of amplitudes considered: 
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from about x = 0.6c, when A = 0, to about 0.34c when A = 0.5 (Figure 
36). 


It is clear that, at this relatively low frequency, the behavior of the 

boundary layer is strongly influenced by the adverse pressure gradients 

present during the first half cycle. The alleviation of these gradients, 

during the second half cycle, has less effect than it does at higher 

frequencies, although the effects are still noticeable. This can be 

seen by comparing the movement of the reversal point for the case: 

U t /c = 4.0, A = 0.5, in Figure 36, with the dashed curve which corre- 
o p 

sponds to the frozen flow: U Q t^/c = 2.0, f^ = 0; the two flows suffer 

the same retardation during the interval 0 < U t/c < 2.0, but different 

— o — 

reta rda t i ons thereafter. 
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1.3 THE NATURE OF UNSTEADY SEPARATION 


Some Important conclusions concerning the nature of unsteady separation 
of the turbulent boundary layer are suggested by the results of the compu- 
tational experiments presented above. 

First: the relaxation of a turbulent boundary layer, which has been 

subjected to an initial retardation strong enough to cause reversal, 
towards a steady state separation, is achieved in a finite time. It was 
suggested, in Reference [1], that the relaxation time is of the same order 
as the time taken for the fluid to be convected off the part of the sur- 
face ahead of the ultimate point of separation. The present results 
support this hypothesis: the calculations for the frozen flow in which 

U Q t^/c = 0.2, f^ = 0, indicate incipient separation at x = 0.27c in a 
time U Q t/c = 2.0, implying that all of the fluid with an average velocity 
exceeding 0.1 5 U q had been replaced by new fluid entering the boundary 
layer after time t^. In the case, just mentioned, the initial retarda- 
tion was abrupt; in cases where the retardation is more gradual the time 
necessary to reach incipient separation is longer. However, the final 
divergence towards singular conditions, as separation is approaches, is 
rapid, with the displacement thickness increasing roughly exponentially 
with time. 

Second: the approach to the separation singularity, in oscillatory flows, 

appears to be a function of both the 1 ong-t i me-average , and the instan- 
taneous, pressure gradients. At low frequencies, the effects of the 
adverse instantaneous pressure gradients seem to dominate, and separation 
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can occur during the first half-cycle if the gradients are sufficiently 

severe. At higher frequencies, the approach to separation appears to 

depend more on the time-average pressure gradients, while the nature of 

the instantaneous gradients may be such as to have a delaying effect on 

separation onset. In one of the cases examined (U t /c = 2.0, A = 0.5), 

o p 

three complete cycles were completed without evidence of imminent separa- 
tion even though the average gradients were large enough to provoke sepa- 
ration after two cycles (as can be seen from the results for U t /c=2.0, 

op’ 

A = 0) . 

The delay of separation onset, observed in the oscillatory cases, where 
the amplitude in large, is evidently associated with the alleviation of 
the pressure gradients during the parts of the cycle when the external 
velocities are increasing with time: instantaneously favorable gradients 

can exist, at sufficiently high amplitudes and frequencies, which offset 
the time-average adverse gradients. On the other hand, the existence of 
favorable instantaneous pressure gradients is not enough to guarantee the 
avoidance of separation -- or even to guarantee the reestablishment of 
forward flow in place of a reversed- f 1 ow region. 

In some cases the reversed- f 1 ow region did vanish during the phases of 
the motion when the pressure gradients were favorable. In others, not 
only did the reversed- f 1 ow region persist during such phases, but the 
extent of the region actually increased. This type of behavior seemed 
to occur when the flow in the boundary layer was already retarded to a 
degree where separation could be expected to occur shortly: the change 

to favorable pressure gradients was then "too late" to save the boundary 
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layer from separating. Figures 37 through 40 contrast, in more detail, 

the cases where recovery was achieved, with those where it was not 

achieved. Figures 37,38 show a typical recovery from reversed flow to 

forward flow, as the pressure gradient changes sign from positive to 

negative. These results, which are for the oscillatory flow: U t /c = 

o p 

2.0, A « 0.5, show the velocity profiles before and after the recovery. 

The boundary- layer thickness increases as a result of the recovery, 

although the displacement thickness decreases (Figure 15). Figures 39, 

40 show a case which fails to recover (U t /c = 4.0, A = 0.5) despite a 

o p 

similar change in the pressure gradient; instead, the reversal point 
moves further upstream, accompanied by an abrupt increase in the boundary- 
layer thickness over the rever sed-flow region. This latter effect is 
reminiscent of the flow patterns, reported in Reference [7], depicting 
the approach to singular conditions, although the onset of a singularity 
(according to the surrogate criterion: 6* = 0.1c) has not yet been 

reached in the solution plotted in Figure 40 . A more detailed examina- 
tion of the results indicates, however, that the solution in Figure 40 
already contains a singularity. Figure 41 shows a plot of the velocity 
at y = 0.1, versus x, for the two sets of results, before and after the 
change in pressure gradient. In the case where recovery is achieved, the 
progression from reversed to forward flow is accompanied by a smooth, 
continuous variation of velocity with x. In contrast, in the case where 
recovery did not take place, a d i scont inu i ty in velocity occurs, and the 
part of the solution where U > 0 loses contact with the part where U < 0. 
The large values of BU/3x, across the discontinuity, give rise to large 
values of 3V/3y, and correspondingly large values of V, which are 
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reflected in the growth of the boundary-layer thickness seen in Figure 
40, but they have not yet been translated into values of 6* exceeding 
0.1c. Nor does the solution imply an excessively large normal pressure 
gradient, at this stage. Two conclusions therefore seem to emerge: 

(1) that the failure of a boundary layer to recover forward flow upon a 
change in the pressure gradient is indicative of an incipient singularity, 
and (2) that the onset of a singularity may not always result in an immediate 
gross thickening of the boundary layer. This second conclusion needs to 
be examined more carefully in future work because singularity onset has 
normally been identified with imminent separation. The present computa- 
tional framework would not be particularly appropriate for such an examina- 
tion because of the coarseness of the mesh in the x-y plane, and because 
there may be some suppression of large temporal gradients (even though 
short time steps were used). Some analytical studies would appear to offer 
more prospect of success. 

Meanwhile, the possibility exists that a sequence of events takes place, 
prior to the final separation of the flow, in which the extent of flow 
reversal increases, and thickening of the boundary layer occurs, regard- 
less of any alleviation of the separat ion-provoki ng pressure gradients. 

In none of the cases discussed here was a retreat from imminent separation 
observed, although it may be demonstrated by further work that the imposi- 
tion of sufficiently strong favorable gradients would bring this about. 

It is interesting to draw a parallel, here, with the hysteresis effect 
observed in dynamic stall: the recovery of the flow on the airfoil is 

not achieved until the incidence is reduced substantially below the level 
reached before stall occurred [10]. The flow around a pitching airfoil 
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is clearly more complicated than the flows discussed in this report; in 
particular, the vortex shed from the leading edge of the airfoil plays an 
important role in the hysteresis process. Nevertheless, a general rule 
may be emerging that, once the separation mechanism is initiated, it is 
necessary to take drastic steps, in terms of a change of external-velocity 
conditions, to restore the boundary layer to an attached condition. Some 
further calculations, on the present lines could usefully be done to pur- 
sue this particular objective. 
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PART II, LAM I NAR FLOW CALCULATIONS 


I 1.1 NATURE OF THE FLOWS CONSIDERED 

The laminar-flow studies were carried out for incompressible, two-dimen- 
sional, time-dependent flow over a flat surface. The viscous portion of 
the flow was taken to be boundary- 1 ayer-1 i ke , insofar as the longitudinal 
extent of the flow was assumed to be large compared with its transverse 
extent, and the potentia 1 -f low velocity was prescribed along the upper edge 
of the integration domain, which ran parallel to the wall. On the other 
hand, no assumption was made about the smallness of second derivatives in 
the longitudinal direction (except at the downstream end of the domain), and 
no restriction was placed upon the variation of static pressure normal to 
the wall: the solution in the interior of the domain corresponds to the 

full unsteady Navier-Stokes equations. Details of the equations, and the 
numerical scheme used to solve them, are given in Appendix B. 

It was assumed that the upstream edge of the plate, upon which the flow 
was developing, was a sharp leading edge. At the downstream end of the 
plate, it was assumed that (regardless of its earlier condition, further 
upstream) the flow had recovered to a form which was consistent with the 
first-order boundary-layer approximations; this assumption was made in order 
to avoid the need for specifying the downstream velocity and vorticity pro- 
files. The region of interest was the portion of the flow over the central 
part of the plate, where reversal was induced by suitable choice of the 
external velocity distributions. These distributions were similar to those 
defined by Equations (1,2,3), above, and the same values of x q /c and x^/c 
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were used. The Reynolds number: U q c/v, was taken to be 10 . 

Results are presented here for a single case, viz. a "frozen flow," of 
the type considered in Part 1, above with U Q t f /c = 2.0, and f^_ = 0.5. 

These results are of interest both in themselves and insofar as they may 
be compared and contrasted with the results of the corresponding turbu- 
lent calculation which has already been described. 

The initial conditions, at t = 0, correspond to steady, constant-pressure 
flow, and were generated by means of a time-relaxation calculation; an 
approximate solution was assumed, and then the calculation was advanced in 
time until the steady-state solution was obtained; At the high Reynolds 
number, for which the calculation was performed, it would be expected that 
the exact solution of the Nav i er-Stokes equations would be close to the 
Blasius solution of the laminar boundary- 1 ayer equations, the only signifi- 
cant difference in boundary conditions lying in the confined domain in 
which the present flow was assumed to be developing; the height of the doma 
was taken to be 0.1c at t = 0. The velocity profiles produced by the time 
relaxation process are compared with the Blasius profiles in Figure 42. 
Calculations were done for three different mesh densities: 

(a) Ax = 0.1c, Ay = 0.01c (11 x 11 node points) 

(b) Ax = 0.05c, Ay = 0.005c (21 x 21 node points) 

(c) Ax = 0.025c, Ay = 0.0025c (41 x 41 node points) 

It will be seen from Figure 42 that the agreement with the Blasius solu- 
tion is good for the finest mesh, tolerable for the medium-density mesh, 
and poor for the coarse mesh, particularly over the upstream half of the 
plate. There are two main reasons for the poor agreement, in this last 
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case: first, the number of node points through the boundary layer is 

too small to give proper resolution of the higher derivatives in the 
y-direction, and second, the wide mesh spacing in the x-direction leads 
to inadequate resolution of the longitudinal derivatives associated with 
the flow near the leading-edge. 

Some comments need to be made here concerning the artificial viscosity, 
introduced by the finite mesh size, and its effect on the quality of the 
solution. The appropriate measures of artificial viscosity are the two 
cell Reynolds numbers: UAx/v and VAy/v, which multiply second deriva- 

tives in the x- and y-d i rect i ons , respectively. If these cell Reynolds 
numbers are everywhere small compared to 2, the solution can be regarded 
as valid; otherwise the effective Reynolds number of the solution is 
reduced, at least in the regions of the flow where those conditions are 
not satisfied. Now, in the present flow, at t = 0, the transverse cell 

Reynolds number: VAy/v is at most about unity, even for the coarse mesh, 

2 2 

and so its effects on the derivatives: 3 /3y , cannot be considered 

serious. The longitudinal cell Reynolds number: UAx/v, is much higher 

than 2, even for the finest mesh, and can scarcely be reduced to 2 without 
placing an enormous number of node points in the x-direction. However, 
it affects only the second derivatives in x which, fortuitously are small 
for the constant-pressure flow; these derivatives are neglected altogether 
in the Bl.asius solution. Consequently, the solution for t = 0 appears 
to be relatively free of the effects of a rt i f i c i a 1 -v i scos i ty . The same 
cannot be said of the subsequent solutions for t > 0, however. There, 
VAy/v would begin to exceed 2, even for the finest mesh, and, also, 
second derivatives in x start to become significant. Thus artificial- 
viscosity becomes a more serious problem as time increases, and this fact 
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needs to be borne in mind during the discussion of the results throughout 
this section. 

On the basis of the comparisons shown in Figure 42, and for reasons of 

computational economy, it was decided to use the medium-density mesh for 

the time-dependent calculations for t > 0. In order to allow for the 

growth of the viscous region, the height of the domain was allowed to 

increase with time: from 0.1c, at t a 0, to 0.2c, at t = 5c/U , after 

o 

which the height was held constant. This scheme proved to be satisfac- 
tory, in the light of the solutions obtained, but it must be remembered 
that the solutions are not independent of the extent of the domain because 
the Navi er-Stokes equations permit pressure gradients to develop normal 
to the wall. The specified rate of increase of height must, therefore, 
be included among the significant boundary conditions imposed on the 
flow. 

The prescribed variation of U, along the upper edge of the domain — 
according to the definition of the frozen flow considered -- represents 
a second boundary condition, and the variation of V represents a third. 

It emerged that the proper variation of V, along the upper edge, had to 
be chosen with some care. I f V is taken to be zero, there, continuity 
demands that the retardation of U, at the boundary, has to be offset by 
an acceleration within the domain. Consequently, the flow near the wall 
feels less retardation than the flow at the upper boundary, and the alle- 
viation of the longitudinal pressure gradients, due to the effects of 
the normal ones, increases with the height of the domain. A retardation 
imposed at a large distance from the wall has little or no effect on the 
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viscous flow, and may not be large enough to produce reversal. In an 
attempt to reduce the extent of the accelerations within the domain, V 
was given a positive value, on the upper boundary, for 0 < x < x q , and 
a negative value for X Q < x < c. Specifically, V was put equal to 
s3U/3x, where s is the height of the domain, which extracted enough mass 
to offset the decrease in U at the boundary. It was not large enough to 
prevent accelerations entirely, within the domain, but experimentation 
showed that such large values of V were needed to achieve this that the 
cell Reynol ds-number problem was seriously aggravated. It suffices to 
say that the performance of the present laminar calculations was not as 
straightforward as it might seem, at first sight. Three seemingly sepa- 
rate problems: the choice of the height of the domain, the choice of 

the V-distribution along the upper boundary, and the art i f i c i a 1 -v i scos i ty 
problem, are actually closely related (computer time is a fourth related 
problem), and a suitable compromise had to be found that attempted to 
sat i sfy all of them. 
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11.2 COMPUTATIONAL EXPERIMENTS 


The calculations for the frozen flow (in which U t Jc = 2.0) were con- 

o T 

tinued up to a time level: U Q t/c “ 8.0. Typical velocity profiles 

through the viscous-flow region are shown in Figures 43 through 45, and 
the longitudinal distributions of wall shear stress and displacement 
thickness are plotted in Figures 46,47- Several of the velocity profiles 
exhibit an overshoot which reflects the acceleration within the doma i n , 
mentioned above, which tends to alleviate the applied retardation. Never- 
theless, in spite of this degree of alleviation, reversal occurs at the 
wa 11 , and the region of reversal extends over most of the plate at the 
highest time level reached. 

Figures 46,47 may be compared with Figures 6,7 which relate to the corre- 
sponding turbulent boundary layer. A number of significant differences 
may be noted between the two sets of results. First, reversal takes 
place before t = t^, in the laminar flow, whereas it did not occur until 
some time later in the turbulent flow. Second, the extent of reversed 
flow is substantially larger in the laminar case. Third, the variation 
of wall shear stress is characteristically different in the two cases: 
it is not only more gradual, in the laminar flow, but it exhibits a double 
minimum at the higher time levels: a double minimum did not occur in the 
present turbulent flows, and indeed has never been observed in any of 
the turbulent calculations performed by us so far. Finally, and perhaps 
most important. The abrupt changes of displacement thickness, associated 
with singularity onset in the turbulent boundary layer, are completely 
absent in the laminar flow. Figure 47 indicates that the displacement 
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thickness is tending towards a stable, steady-state solution; it must be 
remembered that the Navi er-Stokes equations admit solutions corresponding 
to steady-state separation, whereas the first-order boundary- layer equa- 
tions do not. The trend towards an asymptotic solution is emphasized in 
Figures 48,49 which show the variation of the wall shear stress and dis- 
placement thickness with time, for one point on the plate. It is inter- 
esting to note, from Figure 48, that the wall shear stress in the laminar 
flow is of substantially larger magnitude than that in the turbulent case. 
This difference is partly due to the difference in the Reynolds number: 

Zf 

laminar, flat-plate skin friction at R^ « 10 is roughly twice the turbu- 
lent, flat-plate skin friction at R = 10^. The difference is accentuated, 

e 

here, because of the larger extent of the reversed- f low region, in the 

laminar case, and because of the nature of the velocity gradients within 

4 

the region. The displacement thickness in a laminar flow at R^ * 10 is 
about the same as that in a turbulent boundary layer at R^ = 10^, and 
Figure 49 shows that 6“ is of comparable magnitude, in the two cases 
considered here, over the initial time steps. The major difference in 
behavior for larger times, lies in the abrupt increase of 6* with time, 
in the turbulent flow, compared with the progressive increase towards an 
asymptote in the laminar flow. 

It might have been expected that the laminar calculations would exhibit 
some development which, although not of the proportions of a singularity, 
would in some way resemble the steepening gradients characteristic of the 
singularity observed in the turbulent bounda ry- 1 ayer calculations. How- 
ever, the velocity profiles in Figures 43 through 45 show no such devel- 
opment: the abrupt changes in the x-direction — typified by Figure 40, 
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for example -- seem to be completely absent in the present laminar 
results. Inclusion of the second derivatives in x appear to have had the 
effect of inhibiting the formation of any region of unusually large gra- 
dients. The only fear is that art i f i c ia 1 -v i scos i ty effects may have had 
too large an inhibiting effect on those derivatives, and that a more 

accurate solution of the Nav i er-Stokes equations, for a Reynolds number 

4 

of order 10 , might indeed exhibit features somewhat intermediate between 
those of the laminar and turbulent results obtained here. Further 
studies, perhaps in which a large number of node points could be concen- 
trated within a limited region of the flow field, might serve to resolve 
these fears. 

Be that as it may, another area of interest, in the laminar calculations, 
involves the distribution of vorticity throughout the flow, and its 
variation with time. Figures 50 through 52 show contours of constant 
vorticity at three time levels. In the initial phases of the motion, 
the vorticity contours lie essentially parallel to the wall. Subse- 
quently, with increase of time, the contours become distorted, and a 
tongue of h i gh-vor t i c i ty fluid moves out into the body of the flow, 
leaving a region of low or negative vorticity between it and the wall. 
Figure 52 may demonstrate, on a small scale, how vorticity is shed from a 
boundary layer, near separation, and moves away from the body surface. 
There is, of course, no separation "point" as- such; separation is a pro- 
gressive phenomenon, associated with steadily thickening regions of vis- 
cous flow, and correspond i ng modification of the outer inviscid flow. 


31 



CONCLUDING REMARKS 


The calculations reported herein: on the one hand, using the turbulent 

boundary-layer equations, and, on the other, using the laminar Navier- 
Stokes equations, were intended to shed some light on the phenomenon of 
unsteady bounda ry- 1 ayer separation. 

Earlier work had shown that flow reversal, in time-dependent boundary 
layers, is not a singular event, and is not associated with the break- 
down of the boundary-layer approximations. However, it was also demon- 
strated that a singularity could develop some time later: associated with 

steepening streamwise gradients, gross thickening of the boundary layer, 
and the breakdown of those approximations, and that the onset of the 
singularity could be identified with the physical separation of the flow. 
The present calculations of unsteady turbulent boundary layers, which 
included flows relaxing in time and oscillatory flows, confirm the exis- 
tence of a flow regime in which reversal, but not singular conditions, 
prevails. However, they also suggest that singularity onset may not 
always correlate with the immediate gross thickening of the boundary layer 
characteristic of separation. In some of the cases examined, it appeared 
that singularity onset was followed by a phase of only progressive 
boundary- layer thickening; indeed the only unusual gross feature of the 
solution was a lack of responsiveness to changes of the pressure gra- 
dient. In this situation, alleviation of the adverse pressure gradients 
did not bring about a contraction of the reversed- f 1 ow region and the 
rees tabl i shmen t of forward flow over the surface; in some cases the 
extent of the reversed-f low region actually increased. The phenomenon 
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described was most noticeable in oscillatory flows, when reversal was 
found to persist throughout the cycle, despite the presence of instan- 
taneously favorable pressure gradients. 

It is, of course, of dubious validity to draw conclusions from the 
behavior of a solution after a singularity, has developed, and there is 
no direct evidence that the relatively slow increase of boundary- 1 ayer 
thickness, accompanied by the insensitivity to changes of the pressure 
gradient, has any analog in a real flow. On the other hand, there is 
an interesting parallel that can be drawn with the histeresis effect 
observed on airfoils undergoing dynamic stall. 

The laminar-flow calculations, based on the Nav i er-S tokes equations, were 
performed so as to study separation solutions which would not be contami- 
nated by singularity onset and the breakdown of the boundary- 1 ayer equa- 
tions. The solutions obtained certainly had this property, but, sur- 
prisingly, they did not exhibit any of the features associated with the 
early stages of singularity onset: steepening gradients or rapid thicken- 

ing of the viscous region. The calculations were performed for a flow 
relaxing towards a steady-state separation, and the development of the 
solution was almost disappointingly uneventful. The rever sed-f 1 ow region 
was found to increase steadily in extent, the displacement thickness 
increased progressively, approaching an asymptote, and vortical fluid 
moved away from the wall into the interior of the flow. 

The major enigma posed by these results concerns the physical significance 
of the events leading up to singularity onset in the unsteady boundary- 
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layer equations. It was to be expected that the singularity itself was 
simply a feature of the bounda ry- 1 ayer approximations. It was somewhat 
less expected that the formation of a region of locally steep gradients, 
and rapid thickening of the layer, might also prove to be no more than 
a figment of the boundary- 1 ayer imagination. 

The matter is not yet resolved, however, and further work needs to 
address a number of specific questions. First, what is the significance 
of a singularity whose only observable effect is to dissociate the subse- 
quent boundary- layer development from the imposed pressure gradients? 
Second, to what extent was the uneventful development of the laminar 
solutions, obtained here, a result of excessive artificial viscosity, and 
would it genuinely be that uneventful in a high Reynol ds-number viscous 
flow? Perhaps it would be less uneventful in a flow of larger extent 
than the one considered here. Third, at what point, in the introduction 
of higher-order terms into the boundary- 1 ayer equations, does the transi- 
tion from a s i ngu 1 ar i ty-domi na ted separation to a nonsingular one occur, 
and is that transition the same for a laminar flow as it is for a turbu- 
lent flow? 
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APPENDIX A 


CALCULATION METHOD FOR TIME-DEPENDENT 
TURBULENT BOUNDARY LAYERS 


1. Governing Equations 

The velocity components in the x y z - directions are 
expressed in the form U + u, V + v, W + w, respectively, where U,V,W are 
defined as ensemble average velocities, with W 5 0 for a two-dimensional 
flow, and u,v,w are the residual fluctuating components about those 
ensemble averages. The governing equations are similar to those used 
in earlier studies [1,7,9,11], namely, the momentum equation: 


DU + 1 _3jp + JL 

Dt P 3x 3y 


(uv) 


= 0 


the continuity equation: 


(Al) 


+ i V = 0 

3x 3y 


(A2) 


and the empirically-modified turbulent kinetic-energy equation: 


+ 2a, IM TT* * +|r (* 2 33 + r l uv | l/21 * 0 (A3) 


Du v 
Dt 


3U 

3y 


In Equations (Al ,A3) the convective derivative is defined by 


D_ 

Dt 



3x 


v f- 

3y 


(A4) 


In Equation (A3) the empirical convective constant, and the empirical 

functions of position through the boundary layer, and L, are taken to 
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be the same as those used in Reference [1]; accordingly, the same caution- 
ary statements apply to their continued use. Experimental substantiation 
of the appropriateness of Equation (A3) and of the empirical parameters 
in it - or guidance as to how they should be modified or improved - are 
still urgently needed and awaited. The quantity 4> in Equation (A3) is 
the two-dimensional equivalent of the functions 4> x , $ 2 in Reference [12], 
and is defined by 

$ = r{ |uV| + uv ||^-| } , (A5) 

where F is some large number. The inclusion of this term fias no effect 
on the resultant shear stress, but serves to maintain directionality of 
the shear-stress vector according to 

— ~2 ,3 IK 

uv = aj q sgn (— ) , (A6) 

“2 

where q is the resultant mean-square velocity fluctuation. Equation 
(A6) expresses an assumption which is implied in all the two-dimensional 
applications of Townsend's structural similarity hypothesis (see Ref. 

[12]), on which the present model is based. 

2. Solution of the Equations 

The solution of the governing equations follows the approach of 
References [1,7,9]- The flow is divided into an inner and an outer layer, 
with the matching station between them lying at about y = 0.05 . In 

the outer layer the equations are integrated by an explicit, staggered- 
mesh finite-difference scheme, advancing in the positive x-direction. 

The only aspect of the finite-difference scheme, worth mentioning here 
concerns the method of obtaining x-der ? va t i ves . In order to avoid 
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violation of the relevant zones of dependence, derivatives with respect 
to x are formed using two-point backward differences in region of posi- 
tive U, and two-point forward differences in region of reversed flow 
(U < 0) . At stations where the local flow direction is ambiguous: i .e. 

where the sign of U at some node point (x,y,t) is different from that at 
the adjacent points: (x - Ax,y,t) and (x + Ax,y,t), the x-der i vat i ves 

are set equal to zero. This refinement leads to improved smoothness of 
the solution at points of incipient flow reversal. Little loss of accu- 
racy results because U, which is multiplied by the x-derivative in the 
momentum equation, is inevitably close to zero at such points. It should 
be stressed that 3U/3x is not set equal to zero throughout the reversed- 
flow region, as has been done in certain other analyses of flows with 
reversal; such a procedure clearly leads to an invalid solution because 
typical values of U3U/3x are by no means numerically small compared to 
the other terms in the momentum equation. 

Further details of the numerical scheme, in the outer layer, can be 
found in References [1,9,12,13,1**], and the reader is referred thereto. 

In the inner layer, near the surface of the plate, the numerical 
solution is again matched to an approximate solution based on the Law of 
the Wall. The details of this inner solution have been modified, as 
compared with the earlier work, in order to handle the transition from 
positive to negative values of the wall shear stress. In the present 
work, the turbulent kinetic-energy equation, for the inner layer, is 
wr i tten 


3U _ 

3y 


UV I 
Ky 


1/2 


+ A (x ,y , t) } sgn (uv) , 


(A 7) 
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where the dissipation length, L, has been equated, in the usual way, to 
<y. The function A, where 


A = 


1 

2a jUv 


Duv J_ 
Dt uv 



(A8) 


represents the residual imbalance between production and dissipation of 
turbulent kinetic energy which, near the wall, corresponds chiefly to con- 
vective transport. In the inner layer A « |9U/9y|, and is replaced by 

its value, A at the matching point with the ou ter-sol ut ion domain. With 
m 

now independent of y, we have, upon formal integration of Equation (A7) 



o 


uv 

Ky 


1/2 


+ A } sgn 
m 3 


(Uv) dy . 


(A9) 


The integral is evaluated, following Townsend [15], by prescribing a 
linear stress relationship 


- uv = T^/p + ay , (A10) 

where a is independent of y, and typically of the same order as dp/dx. 

It is not difficult to perform the integration analytically; however 
the resulting forms: one, if uv is of the same sign throughout the inner 

layer, and two, if it changes sign, do not lend themselves readily to 
programming for the computer. To avoid these problems, Equation (A9) is 
integrated by a simple iterative numerical scheme. With U and uv known 
at the matching point, an approximate value of x^ is estimated from 
which a is determined; uv and U can then be found, as functions of y, 
for 0 < y < y . The value of x is adjusted, by means of a simple pre- 
dictor-corrector method, until the associated values of U merge smoothly 
into those for the outer domain. 
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3. Integration Domain 


The integration domain extended from x = 0 to x = c , and from 
y = 0 to y = s(x) where s = 1 .255, approximately . The collocation 
points in the y-direction, 20 in number, were distributed, as in Refer- 
ence [13], to give increased density near the wall. The collocation 
points in the x-direction were also distributed nonuniformly, so as to 
give increased density in the center portion of the plate. Specifically, 
the 24 points were distributed according to 

X. = 1.5 t {1 - t + 1 (i-) 2 }, (All) 

j n n 5 n 

with n = 24 and 0 <_ j <_ n . 

4. Boundary Conditions 

For times t £ 0, the flow corresponds to steady, constant-pressure 
flow in the x-direction. The boundary layer is in constant-pressure 
equilibrium, with a thickness, 5, at x = 0, of 0.00444 c. The Reynolds 
number: U q c/v, is taken to be 10^. 

At x = 0, the velocity and shear-stress profiles are maintained, 
for all time, in the same form as at t = 0; i.e., in steady constant- 
pressure equilibrium with an external velocity of and the Reynolds- 
number conditions specified above. 


39 



APPENDIX B 


SOLUTION METHOD FOR THE T I ME -DEPENDENT NAV I ER-STOKES EQUATIONS 
1» Governing Equations 

The Navier Stokes equations, for incompressible, time-dependent 
flow, can be written in vort ici ty-transport form as 


DC _ K 
Dt at 




(Bl) 


in which the single non-zero component of vort i city is where 


3V _ 3 U_ 

3x 3y * 


and 


2 2 
v 2 = iL + iL 
-2 .2 
ax ay 


From Equation (B2) , together with the continuity equation: 


(B2) 


(B3) 


au + av 

ax ay 


o , 


(BA) 


a Poisson equation can be derived for the velocity components, in terms 
of derivatives of the vorticity: 


V 2 U = 


is. 

sy 


(B5) 


v 2 v 


K 

ax ' 


(B6) 


Equations (B1,B5,B6) form a parabolic set, and represent the governing 
equations, in three unknowns: C,U,V, which are integrated in a three- 



dimensional domain consisting of two space dimensions and time. The 
method of solution is closely related to that of References [16,17] for 
solving the parabol icized Nav i er-Stokes equations in three-dimensional 
steady flow. 

2. Solution of the Equations 

The governing equations are integrated by means of an implicit, alterna- 
ting-direction (ADI) scheme, advancing in the positive time direction. 

A rectangular mesh: x = mAx, y = nAy, is erected on the plate, permitting 

discretization of the derivatives at each time level, l: 


where F 


3F _ 

F U) . F U) 

1 , n m , n 

if U < 0, 

(B7) 

3x 

Ax 


F U) . F U) 

m, n m- 1 , n 

, if U _> 0; 

(B8) 


Ax 

3F _ 

F U) . F (i) 

m,n+l m,n 

if V < 0, 

(B9) 

3y 

Ay 


pU) _ p(A) 
m , n m , n- 1 

if \i > n • 

(BIO) 


Ay 


1 

Ll. 

CM 

CO 

F U) _ 2F (i) + 

m+l ,n m,n 

F U) 

m-1 ,n 

(Bll) 

3x 2 

(Ax) 2 

f 

3 2 F . 

. - 2F^ + 
m,n+l m,n 

F U) 

m,n-l 

( B 1 2 ) 

a? 

*< 

ro 

1 

(Ay) 2 

9 

3F _ 

F U) . F U-1) 

m , n m , n 


(B13) 

at 

At 


is any variable. 

Any one of the 

governing equations 

can then be 
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written in difference form as 


AF i+! + BF ^+? + CF ^ + DF ^! + EF ^ 

m+I ,n m,n+i m,n m-l,n m,n-l 


R 


(Bl4) 


where A through E are coefficients which, in general, depend on the solu- 
tion, but which are regarded as known at each iteration level reflecting 
the customary linearization procedure. R, in Equation (B14) involves 
the solution at time level £-1, which has already been calculated. The 
field, m,n, is scanned alternately in the m-, and n-di rect ions , convert- 
ing Equation (Bl4) into the successive forms: 


AF 


BF 


U) 

m+1 ,n 

+ CF^ + DF^J = R, 

m,n m-1 ,n 1 

(B 1 5 ) 

U) 

m,n+l 

+ CF^ + EF^ , = R, 
m,n m,n-l 2 

(B 1 6 ) 


where Rj and have absorbed the passive terms on the left-hand side 
of Equation ( B 1 4) . Equations (B15»B16) are solved by the extended 
Choleski method [18] , for which efficient solution algorithms are 
ava i lab 1 e. 


The procedure during any iteration cycle is to solve the vorticity equa- 
tion to provide updated values of £, and then, in a second step, to solve 
the two Poisson equations to provide updated values of the velocities 
U,V. The Poisson equations were solved simultaneously, by regarding F 
and R, in equation (B 1 4) as two-d imens iona 1 vectors, and A through E as 
square matr i ces . 

3. Boundary Conditions 

It is assumed that the upper surface of the integration domain (y = s) 
is outside the viscous region, so that £ = 0. The velocity components 
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are assumed to be prescribed. At the upstream boundary (x = 0) , the 
vorticity is assumed to be zero except at the wall where a delta func- 
tion is imposed, corresponding to the no-^lip condition. At the down- 
stream boundary (x « c) , the vorticity is assumed to be zero outside 
the boundary layer, while conditions inside the boundary layer are 

assumed to conform to the first-order boundary- 1 ayer approximations: 

2 2 

9 /9x = 0. At the wall: y = 0, the boundary conditions are U = V = 0; 

the vorticity at the wall is determined as part of the solution, and 
is proportional to the wall shear stress. 
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Figure 4 Wall Shear Stress Distributions for Increasing Time; Frozen Flow, Ugt^/c = 2.0, = 0. 
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Figure 6 Wall Shear Stress Distributions for Increasing Time; 
Frozen Flow, t £ = 2.0, = 0.5. 
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Figure 16 Wall Shear Stress Distributions for Increasing Time; 
Oscillatory Flow, U Q tp/c = 4 . 0 , A = 0.5. 
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Figure 22 Variation of Displacement Thickness with Time; Oscillatory Flow, 
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U t_/c = 2.8. 
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Figure 25 Wall Shear Stress Distributions for Various Amplitudes; 
Oscillatory Flow, U t^/c = 2.0, U Q t/c =4.0. 











Figure 30 Variation of Displacement Thickness with Time; 

Oscillatory Flow, U t /c = 2.0, x = 0.65c. 
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Figure 34 Wall Shear Stress Distributions for Various 

Amplitudes; Oscillatory Flow, U t /c = 4.0, U t/c = 4.0. 
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Figure 39 


Velocity Profiles at U Q t/c = 



2.0; Oscillatory Flow, U t /c = 4.0, A = 0.5. 
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Figure 40 


Velocity Profiles at U^t/c = 3.0; Oscillatory Flow, 


U Q tp/c = 4.0, A = 0.5. 
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Figure 49 Variation of Displacement Thickness with Time; 

Frozen Flow, l^t^/c = 2.0, = 0.5, x = 0.65c. 






